**** Honey experiment****

/*

Kids first

*/

if 0{
* 0. load data for kids

use "${dataanalysis}\KidsReplication.dta", clear

encode id, gen(idcode)
xtset idcode post
bysort idcode (post): gen worked = dem_work[1]
* 1. introductory things and globals

gen lnprice = ln(wtp_honey)
la  var wtp_honey "Price offered honey"

global controlskids dem_gender dem_age  timesincetreatment spillsample simpleindex
set scheme tab1

quietly{ 
	foreach var in $controlskids wtp_notebook_yesno {
	sum `var'
	g `var'_m = r(mean)
	
}
}

fgen lnhoney = ln(wtp_honey)

* 2. Graphs

preserve

collapse (mean) meanpurchased = wtp_honey_yesno (sd) sdpurchased = wtp_honey_yesno (count)  n= wtp_honey_yesno, by(wtp_honey treatindividual dem_gender)
generate hi = meanpurchased + invttail(n-1,0.05)*(sdpurchased / sqrt(n))
generate lo = meanpurchased - invttail(n-1,0.05)*(sdpurchased / sqrt(n))

twoway (scatter meanpurchased wtp_honey if treatindividual ==1 & dem_gender == 0) ///
         (scatter meanpurchased wtp_honey if treatindividual == 0 & dem_gender == 0, lpattern(dash)) ///
		 (rcap hi lo wtp_honey if treatindividual == 1 & dem_gender == 0, color(gs13)), ///	
		 xtitle(Offer price) ytitle(Proportion purchased) legend (position(6) order(1 "Treated" 2 "Control" 3 "90% Confidence Intervals")) graphregion(color(white)) bgcolor(white) ylabel(0(.25)1) xlabel(90(10)130)
		 
		 
graph export "${figures}\honeyboys_simple.eps", replace
	 
twoway (scatter meanpurchased wtp_honey if treatindividual ==1 & dem_gender == 1) ///
         (scatter meanpurchased wtp_honey if treatindividual == 0 & dem_gender == 1, lpattern(dash)) ///
		 (rcap hi lo wtp_honey if treatindividual == 1 & dem_gender == 1, color(gs13)), ///	
		 xtitle(Offer price) ytitle(Proportion purchased) legend (position(6) order(1 "Treated" 2 "Control" 3 "90% Confidence Intervals")) graphregion(color(white)) bgcolor(white) ylabel(0(.25)1) xlabel(90(10)130)
		 
		 
graph export "${figures}\honeygirls_simple.eps", replace
		 		 
 
restore


* 3. Estimation for median wtp

* sum stats

bysort treatindiv dem_gender: sum wtp_honey_yesno

* no controls

probit wtp_honey_yesno lnhoney spillsample if treatindividual == 0
nlcom (WTP:exp(- (_b[_cons])/_b[lnhoney])),noheader post

probit wtp_honey_yesno lnhoney if treatindividual == 1
nlcom (WTP:exp(- (_b[_cons])/_b[lnhoney])),noheader post


* wtp for mean individual in control
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 0
eststo wtpc: nlcom (WTP:exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment]+dem_gender_m*_b[dem_gender] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bc = r(b)
matrix sc = r(V)

probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 1
* wtp for mean individual in treatment group: 
eststo wtpt: nlcom (WTP: exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment]+dem_gender_m*_b[dem_gender] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bt = r(b)
matrix st = r(V)

matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2)) // Not sure this is the right test -- perhaps not appropriate for medians
matrix list stat

* wtp for mean individual in control -- boys
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 0 & dem_gender == 0
eststo wtpcb: nlcom (WTP:exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bc = r(b)
matrix sc = r(V)

probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 1 & dem_gender == 0
* wtp for mean individual in treatment group: -- boys
eststo wtptb: nlcom (WTP: exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bt = r(b)
matrix st = r(V)
matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2)) 

matrix list stat

probit wtp_honey_yesno i.treatindividual##c.lnhoney wtp_notebook_yesno $controlskids if dem_gender == 0
contrast treatindividual treatindividual#c.lnhoney , overall


* wtp for mean individual in control -- girls
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 0 & dem_gender == 1
eststo wtpcg: nlcom (WTP:(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post
matrix bc = r(b)
matrix sc = r(V)

probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controlskids if treatindividual == 1 & dem_gender == 1
* wtp for mean individual in treatment group: -- girls
eststo wtptg: nlcom (WTP:exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bt = r(b)
matrix st = r(V)
matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2))

matrix list stat
probit wtp_honey_yesno i.treatindividual##c.lnhoney wtp_notebook_yesno $controlskids if dem_gender == 0
contrast treatindividual treatindividual#c.lnhoney , overall

probit wtp_honey_yesno i.treatindividual##c.lnhoney wtp_notebook_yesno $controlskids if dem_gender == 1
contrast treatindividual treatindividual#c.lnhoney, overall


esttab wtpc wtpt  using "${tables}\wtpkids.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  replace  fragment collabels(none) ///
	  posthead( \emph{Combined sample} & \\)
	  
esttab wtpcb wtptb  using "${tables}\wtpkids.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  append  fragment collabels(none) ///
	  prehead( ) ///
		posthead(\midrule \emph{Boys} \\ ) 

esttab wtpcg wtptg  using "${tables}\wtpkids.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  append  fragment collabels(none)  ///
	  prehead( ) ///
		posthead(\midrule \emph{Girls} \\ ) 


* 4.  Estimation for more detailed graphs


graph bar wtp_honey_yesno, over(wtp_honey) over(treatindividual) ytitle("Proportion purchased honey")
graph export "${figures}\honeykids.eps", replace


probit wtp_honey_yesno lnhoney
eststo col1: margins,  dydx(*) post
probit wtp_honey_yesno lnhoney i.treatindividual
eststo col2: margins,  dydx(*) post
probit wtp_honey_yesno i.treatindividual##c.lnhoney
eststo col3: margins,  dydx(*) post
probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno  $controlskids
eststo col4: margins,  dydx(*) post



probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno o.wtp_honey  $controlskids
f_able lnhoney, auto
margins treatindividual, at(wtp_honey=(90(10)130)) atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeycontrols.eps", replace

probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno o.wtp_honey  $controlskids if dem_gender == 1
f_able lnhoney, auto
margins treatindividual, at(wtp_honey=(90(10)130)) atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeygirls.eps", replace

probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno o.wtp_honey  $controlskids if dem_gender == 0
f_able lnhoney, auto
margins treatindividual, at(wtp_honey=(90(10)130))  atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeyboys.eps", replace

 }

/*

Now parents

*/

if 1 {
* 0. load data for adults

use "${dataanalysis}\ParentsReplication.dta", clear

* 1. introductory things and globals

la  var wtp_honey "Price offered honey"
replace dem_age = 40 if dem_age == -999
drop if spillsample == 1

global controls dem_gender dem_age simpleindex timesincetreatment spillsample 
set scheme tab1

quietly{ 
	foreach var in $controls wtp_notebook_yesno  {
	sum `var'
	g `var'_m = r(mean)
	
}
}


fgen lnhoney =  ln(wtp_honey)

* 2. Graphs

graph bar wtp_honey_yesno, over(wtp_honey) over(treatindividual) ytitle("Proportion purchased honey")
graph export "${figures}\honeyadult.eps", replace


* 3. Estimation median WTP

* sum stats

bysort treatindiv dem_gender: sum wtp_honey_yesno

* linear probability model

reg wtp_honey_yesno 1.treatindividual wtp_honey  control50 control75
reg wtp_honey_yesno 1.treatindividual i.wtp_honey control50 control75


reg wtp_honey_yesno 1.treatindividual##i.dem_gender wtp_honey control50 control75
reg wtp_honey_yesno 1.treatindividual##i.dem_gender i.wtp_honey control50 control75



reg wtp_honey_yesno 1.treatindividual wtp_honey wtp_notebook_yesno $controls 
reg wtp_honey_yesno 1.treatindividual i.wtp_honey wtp_notebook_yesno $controls



reg wtp_honey_yesno 1.treatindividual#i.dem_gender wtp_honey wtp_notebook_yesno $controls
reg wtp_honey_yesno 1.treatindividual#i.dem_gender i.wtp_honey wtp_notebook_yesno $controls


* wtp for mean individual in control
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 0
eststo wtpc: nlcom (WTP: exp(-(_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment]+dem_gender_m*_b[dem_gender] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bc = r(b)
matrix sc = r(V)


probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 1
* wtp for mean individual in treatment group: 
eststo wtpt: nlcom (WTP: exp(-(_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment]+dem_gender_m*_b[dem_gender] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bt = r(b)
matrix st = r(V)

matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2))
matrix list stat

* wtp for mean individual in control -- men
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 0 & dem_gender == 0
eststo wtpcb: nlcom (WTP: exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bc = r(b)
matrix sc = r(V)

probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 1 & dem_gender == 0
* wtp for mean individual in treatment group -- men
eststo wtptb: nlcom (WTP:exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post


matrix bt = r(b)
matrix st = r(V)
matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2))


matrix list stat

* wtp for mean individual in control -- women
probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 0 & dem_gender == 1
eststo wtpcg: nlcom (WTP: exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post
matrix bc = r(b)
matrix sc = r(V)

probit wtp_honey_yesno lnhoney wtp_notebook_yesno $controls if treatindividual == 1 & dem_gender == 1
* wtp for mean individual in treatment group: -- women
eststo wtptg: nlcom (WTP: exp(- (_b[_cons]+simpleindex_m*_b[simpleindex]+timesincetreatment_m*_b[timesincetreatment] + ///
dem_age_m*_b[dem_age]+ spillsample_m*_b[spillsample]+wtp_notebook_yesno_m*_b[wtp_notebook_yesno])/_b[lnhoney])),noheader post

matrix bt = r(b)
matrix st = r(V)
matrix stat = (bt[1,1]-bc[1,1])/(sqrt(sc[1,1]^2+st[1,1]^2))

matrix list stat




esttab wtpc wtpt  using "${tables}\wtpadults.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  replace  fragment collabels(none) ///
	  posthead( \emph{Combined sample} & \\)
	  
esttab wtpcb wtptb  using "${tables}\wtpadults.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  append  fragment collabels(none) ///
		posthead(\midrule \emph{Men} \\ ) 

esttab wtpcg wtptg  using "${tables}\wtpadults.tex", style(tex) cells(b(star fmt(3)) se(par)) ///
		starlevels(* 0.10 ** 0.05 *** 0.01) label mlabels("" "" ""  "" "") ///
	  nonumbers  append  fragment collabels(none)  ///
		posthead(\midrule \emph{Women} \\ ) 


* 4. Estimation for graphs


probit wtp_honey_yesno lnhoney
eststo col1: margins,  dydx(*) post
probit wtp_honey_yesno lnhoney i.treatindividual
eststo col2: margins,  dydx(*) post
probit wtp_honey_yesno i.treatindividual##c.lnhoney
eststo col3: margins,  dydx(*) post
probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno  $controls
eststo col4: margins,  dydx(*) post

probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno  o.wtp_honey $controls
f_able lnhoney, auto

margins treatindividual, at(wtp_honey=(90(10)130)) atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeycontrolsadult.eps", replace


probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno o.wtp_honey $controls if dem_gender == 1
f_able lnhoney, auto
margins treatindividual, at(wtp_honey=(90(10)130)) atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeywomen.eps", replace

probit wtp_honey_yesno c.lnhoney##i.treatindividual  wtp_notebook_yesno o.wtp_honey  $controls if dem_gender == 0
f_able lnhoney, auto
margins treatindividual, at(wtp_honey=(90(10)130)) atmean noestimcheck
marginsplot, ytitle(Pr(Purchase)) ciopts(lpattern(dash)) plot1opts(lcolor(orange*.8)  mcolor(orange*.8)) ci1opts(lcolor(orange*.8) ) title("Predictive margins by treatment status")
graph export "${figures}\honeymen.eps", replace

}
